The present document contains details about the data-analysis for the paper entitled: “The assessment of presence and performance in an AR environment for motor imitation learning: a case-study on violinists.” Authors: Adriaan Campo, Aleksandra Michałko, Bavo Van Kerrebroeck, Boris Stajic, Maja Pokric, Marc Leman. Basically, the paper tests a violin playback system. The violinist’s task is to play in synchrony with the principal violinist of an orchestra (represented as an avatar). The playback system provides the audio, in addition to a 2D or 3D avatar via a Hololens. The focus is mainly on the effect due to the experimental conditions of a 2D or 3D playback.
Figure 1 shows two statistical workflows. The upper workflow on top
is about the biomechanical metrics from motion capture (here called:
procustus and sparc, instead of PD and dSI as in the paper). Each metric
represents differences between the violinist and the avatar, summarized
over time. These values are used as response in regression models,
whereas conditions, participants and trial are used as predictors.
Model_1 and model_2 are similar models except that in model_2
difficulty is added in interaction with condition. The
lower workflow is about the behavioral metrics of the questionnaires. We
have 3 presence questionnaires. Model_3 and model_4 are similar models
except that in model_4 procustus is added in interaction
with condition.
In the upper workflow, we start with a comparison of model_1 and
model_2 to test whether difficulty should be added to the
model. We then perform a more detailed diagnostics of the best model, as
well as a contrast analysis of condition and trials. In the lower
workflow, we start with a comparison of model_3 and model_4 to test
whether procustus should be added to the model. We then
proceed with a more detailed diagnostics of the best model, as well as a
contrast analysis of conditions.
The workflows hold a scheme for testing the work hypothesis, as shown in figure 2:
Hypothesis 1. Students will show better violin performance in the 3D condition compared to the 2D condition:
1.1. Similarity between virtual teachers’ and students’ bow movement is higher .
1.2. Movement smoothness is higher
Hypothesis 2. Learning effectiveness of violin performance will be higher in the 3D condition compared to the 2D condition.
Hypothesis 3. The 3D condition will induce a higher level of presence compared to the 2D condition:
3.1. “physical presence” will be higher
3.2. “social presence” will be higher
Hypothesis 4. The level of presence in AR influences students’ violin performance.
In hypothesis 1, we first compare model_1 and model_2 and focus on
the contrast analysis of the best model. This analysis will tell us the
differences between conditions. In hypothesis 2, we expand our contrast
analysis by comparing conditions and trials. In hypothesis 3, we do a
contrast analysis of model_3.
In hypothesis 4 we test whether the procustus metric might be a relevant
predictor variable.
A retrospective power analysis of the setup shown in Figure 3 reveals that 11 participants and 4 repeated measures has sufficient power when Cohen’s D is > 0.6. To show that, we used a hierarchical statistical model with 2 conditions, and with participants and trials modeled as groups for which we assumed a standard deviation of 0.25 and 0.1, respectively. We tested Cohen’s D from 0 to 0.6, and a number of participants from 5 to 30. For each dot in the graph below (e.g. D = 0.5, number of participants = 20) we simulated 500 models and tested whether the contrast of conditions has a probability mass of >= .95. The proportion of yes (versus no) is presented as power (in %). The results show that 11 participants have enough power (> 80%) when Cohen’s D is > 0.6, meaning that there is <20% to miss an effect when there is in fact an effect (= false negative). In our models, we observe that the calibrated models have D > 0.8 (for the models: see below).
The first box in the workflow shown in Figure 1 is about data. We use one dataset for all models. The data variables are listed below:
load(file = "Data.RData")
str(Data)
## 'data.frame': 88 obs. of 16 variables:
## $ participant : Factor w/ 11 levels "1","2","3","4",..: 1 1 1 1 1 1 1 1 10 10 ...
## $ condition : Factor w/ 2 levels "1","2": 1 1 1 1 2 2 2 2 1 1 ...
## $ trial : Factor w/ 4 levels "1","2","3","4": 1 2 3 4 1 2 3 4 1 2 ...
## $ response1M_procustus : num 0.5448 0.0614 0.1508 -0.076 -0.0535 ...
## $ response0M_procustus : num 1.319 0.835 0.925 0.698 0.43 ...
## $ response1M_sparc : num -0.0896 -0.1224 -0.1035 -0.048 0.1661 ...
## $ response0M_sparc : num 11.2 11.1 11.1 11.2 11.1 ...
## $ WPQ : num 4.07 3.93 4.04 4.11 5 ...
## $ MPQS : num 2.8 3.2 3 2 4 5.2 4.6 4.8 4.8 4 ...
## $ MPQP : num 3.4 3 2.8 2.8 4.2 6 5.2 6 3.8 4.4 ...
## $ Difficulty : num [1:88, 1] -0.203 -0.993 -1.586 -1.191 -0.5 ...
## $ log_response0M_procustus: num 0.2766 -0.18 -0.0784 -0.3598 -0.8436 ...
## $ log_response1M_procustus: num -0.607 -2.79 -1.892 NA NA ...
## $ log_response0M_sparc : num 2.41 2.41 2.41 2.42 2.4 ...
## $ log_response1M_sparc : num NA NA NA NA -1.8 ...
## $ difficulty_s : num [1:88, 1] -0.203 -0.993 -1.586 -1.191 -0.5 ...
head(Data)
## participant condition trial response1M_procustus response0M_procustus
## 1 1 1 1 0.54477643 1.3185912
## 2 1 1 2 0.06143131 0.8352461
## 3 1 1 3 0.15080784 0.9246226
## 4 1 1 4 -0.07599130 0.6978235
## 5 1 2 1 -0.05350941 0.4301667
## 6 1 2 2 -0.07670672 0.4069694
## response1M_sparc response0M_sparc WPQ MPQS MPQP Difficulty
## 1 -0.08964702 11.15673 4.074074 2.8 3.4 -0.2033571
## 2 -0.12239219 11.12399 3.925926 3.2 3.0 -0.9932761
## 3 -0.10351349 11.14287 4.037037 3.0 2.8 -1.5857154
## 4 -0.04796922 11.19841 4.111111 2.0 2.8 -1.1907559
## 5 0.16611262 11.06114 5.000000 4.0 4.2 -0.4995767
## 6 0.16967439 11.06470 5.518519 5.2 6.0 0.6359318
## log_response0M_procustus log_response1M_procustus log_response0M_sparc
## 1 0.27656388 -0.6073798 2.412043
## 2 -0.18002891 -2.7898356 2.409104
## 3 -0.07836964 -1.8917489 2.410800
## 4 -0.35978914 NA 2.415772
## 5 -0.84358250 NA 2.403438
## 6 -0.89901733 NA 2.403760
## log_response1M_sparc difficulty_s
## 1 NA -0.2033571
## 2 NA -0.9932761
## 3 NA -1.5857154
## 4 NA -1.1907559
## 5 -1.795089 -0.4995767
## 6 -1.773874 0.6359318
summary(Data)
## participant condition trial response1M_procustus response0M_procustus
## 1 : 8 1:44 1:22 Min. :-0.279936 Min. :0.3399
## 2 : 8 2:44 2:23 1st Qu.:-0.088432 1st Qu.:0.5096
## 3 : 8 3:21 Median : 0.008438 Median :0.6013
## 4 : 8 4:22 Mean : 0.056997 Mean :0.6749
## 5 : 8 3rd Qu.: 0.111951 3rd Qu.:0.8061
## 6 : 8 Max. : 1.247398 Max. :2.0212
## (Other):40
## response1M_sparc response0M_sparc WPQ MPQS
## Min. :-0.640874 Min. :10.25 Min. :3.852 Min. :1.000
## 1st Qu.:-0.085935 1st Qu.:11.00 1st Qu.:4.444 1st Qu.:2.600
## Median : 0.005891 Median :11.17 Median :4.907 Median :3.200
## Mean :-0.018900 Mean :11.10 Mean :4.891 Mean :3.151
## 3rd Qu.: 0.090802 3rd Qu.:11.28 3rd Qu.:5.259 3rd Qu.:3.800
## Max. : 0.309727 Max. :11.56 Max. :6.222 Max. :5.200
## NA's :2 NA's :2
## MPQP Difficulty.V1 log_response0M_procustus
## Min. :2.000 Min. :-2.0300448 Min. :-1.0792
## 1st Qu.:3.400 1st Qu.:-0.7217415 1st Qu.:-0.6742
## Median :4.200 Median :-0.0799323 Median :-0.5086
## Mean :3.991 Mean : 0.0000000 Mean :-0.4594
## 3rd Qu.:4.600 3rd Qu.: 0.7470142 3rd Qu.:-0.2156
## Max. :6.000 Max. : 2.2651397 Max. : 0.7037
## NA's :2
## log_response1M_procustus log_response0M_sparc log_response1M_sparc
## Min. :-6.2271 Min. :2.328 Min. :-5.680
## 1st Qu.:-3.3151 1st Qu.:2.398 1st Qu.:-3.371
## Median :-2.6006 Median :2.413 Median :-2.401
## Mean :-2.4148 Mean :2.407 Mean :-2.657
## 3rd Qu.:-1.5421 3rd Qu.:2.423 3rd Qu.:-1.777
## Max. : 0.2211 Max. :2.447 Max. :-1.172
## NA's :41 NA's :43
## difficulty_s.V1
## Min. :-2.0300448
## 1st Qu.:-0.7217415
## Median :-0.0799323
## Mean : 0.0000000
## 3rd Qu.: 0.7470142
## Max. : 2.2651397
##
It is of interest to show the differences between non-calibrated and calibrated metrics, all obtained by taking the median of the data points (summarizing time).
load(file = "Data.RData")
library(ggplot2)
library(patchwork)
p1 <- ggplot(Data) +
geom_point(aes(x=participant,y=response0M_procustus, color = factor(condition)), position = position_jitterdodge()) +
ggtitle("procustus")
p2 <- ggplot(Data) +
geom_point(aes(x=participant,y=response1M_procustus, color = factor(condition)), position = position_jitterdodge()) +
ggtitle("procustus-calibrated")
p3 <- ggplot(Data) +
geom_point(aes(x=participant,y=response0M_sparc, color = factor(condition)), position = position_jitterdodge()) +
ggtitle("sparc")
p4 <- ggplot(Data) +
geom_point(aes(x=participant,y=response1M_sparc, color = factor(condition)), position = position_jitterdodge())+
ggtitle("sparc-calibrated")
p1 / p2
p3 / p4
Here we show the answers to 3 presence questionnaires and 1 question about perceived difficulty.
load(file = "Data.RData")
library(ggplot2)
library(patchwork)
p5 <- ggplot(Data) +
geom_point(aes(x=participant,y=WPQ, color = factor(condition)), position = position_jitterdodge())+
ggtitle("Witmer Presence questionnaire")
p6 <- ggplot(Data) +
geom_point(aes(x=participant,y=MPQS, color = factor(condition)), position = position_jitterdodge())+
ggtitle("social Makransky multimodal presence scale")
p7 <- ggplot(Data) +
geom_point(aes(x=participant,y=MPQP, color = factor(condition)), position = position_jitterdodge())+
ggtitle("physical Makransky multimodal presence scale")
p8 <- ggplot(Data) +
geom_point(aes(x=participant,y=difficulty_s, color = factor(condition)), position = position_jitterdodge())+
ggtitle("Perceived difficulty")
p5 / p6
## Warning: Removed 2 rows containing missing values (`geom_point()`).
## Removed 2 rows containing missing values (`geom_point()`).
p7 / p8
## Warning: Removed 2 rows containing missing values (`geom_point()`).
# p8 <- ggplot(Data) +
# geom_point(aes(x=participant,y=responseTSM_procustus, color = condition), position = position_jitterdodge())+
# ggtitle("procustus-calibrated via smooth regression")
# p8
The original data of the metrics procustus and sparc are time series.
The second box in the workflow of Figure 1 is about regression
modelling. We tested several models but we ended up with four basic
models, model_1 and model_2 for the metric workflow, and model_3 and
model_4 for the questionnaire workflow. The syntax of the models (here
in R package brms format) is very similar.
Model_1: response ~ 0 + condition + (1 | condition:participant + condition:trial) )
Model_2: response ~ 0 + condition*difficulty + (1 + difficulty | condition:participant + condition:trial) )
Model_3: response ~ 0 + condition + (1 | condition:participant + condition:trial) )
Model_4: response ~ 0 + condition*procustus + (1 + procustus | condition:participant + condition:trial) )
where:
response is either the procustus of sparc metrics
(giving us 2 different models for model_1),
condition is either a 2D or 3D rendering of the
visual scene,
difficulty is the perceived difficulty of the
task,
procustus is the procustus metric of the
task,
trial is the participant’s session.
In model_1 and model_2 a skew_normal link function is
used, in model_3 and model_4 a gaussian link function is
used. Note further that participant and trial
are exchangeable variables. The advantage of the mixed model is that
these variables can be modelled as instances of distributions at a
higher hierarchical level. Accordingly, each participant, being
exchangeable, is drawn from a normal distribution whose sd is estimated
by the model. Same for trial. This modelling approach prevents
overfitting by shrinking the instances of the group-level variables
participant and trial towards the means of the
respective group-level. Since condition has only two
levels, we keep it as population variable. Group-level effects of
trial are used later in a contrast analysis. Another way of
looking at this regression is that it captures variability that is
related to participant and trial, leaving a
more “pure” variability of interest to condition.
We run the models on a 48 dual core machine (at Ghent University,
IPEM), using the R package brms. We take 5000 warmups and
40000 iterations, with an adapt_delta = 0.995 and max_treedepth = 12, 4
chains, and 24 threads. The large amount of iterations was needed in
view of a stable Bayes factor test in the R package
parameters.
We then proceed with the analysis in 3 parts (Figure 1).
difficulty, model_2 = with difficulty) using
the Bayes-factor test (using bayes_factor()). Running
ahead, we found that none of the model_2 turn are any better than
model_1.pp_check() for a global
retrodiction check and mcmc_plot() for an overview of the
posterior distributions of parameters, we also run a
bayes_R2() to get an estimate of the variances, and
parameters() in order to get a summary of the model.integer
(rather than factor) but we thought that a factor approach was more
appropriate given the fact that order was relevant, instead of the exact
time between the sessions. We report contrast testing both as table and
as plot.Part 1 of this analysis is related to the procustus and sparc metrics and hypothesis 1 and 2.
We tested the models for procustus and sparc and report here the log of the Bayes factor.
load(file = "Results/BF_model_comparison_procustus0M.RData")
load(file = "Results/BF_model_comparison_sparc0M.RData")
print(paste("Bayes factor in favor of model_1 over model_2 (procutus0M): ", BF_model_comparison_procustus0M$bf ))
## [1] "Bayes factor in favor of model_1 over model_2 (procutus0M): 6.78245557353672"
print(paste("Bayes factor in favor of model_1 over model_2 (sparc0M): ", BF_model_comparison_sparc0M$bf ))
## [1] "Bayes factor in favor of model_1 over model_2 (sparc0M): 278.982285441945"
load(file = "Results/BF_model_comparison_procustus1M.RData")
load(file = "Results/BF_model_comparison_sparc1M.RData")
print(paste("Bayes factor in favor of model_1 over model_2 (procustus1M): ", BF_model_comparison_procustus1M$bf ))
## [1] "Bayes factor in favor of model_1 over model_2 (procustus1M): 1.44338762610906"
print(paste("Bayes factor in favor of model_1 over model_2 (sparc1M): ", BF_model_comparison_sparc1M$bf ))
## [1] "Bayes factor in favor of model_1 over model_2 (sparc1M): 1870.12670438265"
We conclude that there is strong evidence for model_1 (i.e., without
difficulty).
We show the diagnostics both for the procustus and sparc model_1:
the model_1 formula,
the Bayes_R2 analysis,
the model parameters,
the posterior prediction check (pp_check) next to the plot of
fixed parameters (i.e. condition)
load(file = "Results/post_analysis_procustus0M.RData")
load(file = "Results/post_analysis_sparc0M.RData")
post_analysis_procustus0M
## [[1]]
## response0M_procustus ~ 0 + condition + (1 | condition:participant) + (1 | condition:trial)
##
## [[2]]
##
## [[3]]
## R2 SD CI CI_low CI_high CI_method Component Effectsize
## 1 0.55 0.07 0.95 0.4 0.68 HDI conditional Bayesian R-squared
## 2 0.06 0.07 0.95 0.0 0.21 HDI marginal Bayesian R-squared
##
## [[4]]
## Parameter Effects Component Mean CI CI_low
## 1 b_condition1 fixed conditional 0.74 0.95 0.60
## 2 b_condition2 fixed conditional 0.61 0.95 0.47
## 3 sd_condition:participant__Intercept random conditional 0.18 0.95 0.12
## 4 sd_condition:trial__Intercept random conditional 0.08 0.95 0.02
## 5 sigma fixed sigma 0.15 0.95 0.12
## CI_high pd log_BF Rhat ESS Group
## 1 0.87 1 15.50 1.00 963.21
## 2 0.74 1 10.43 1.00 764.86
## 3 0.26 1 8.72 1.00 847.66 condition:participant
## 4 0.17 1 -0.74 1.01 716.15 condition:trial
## 5 0.18 1 26.52 1.00 1196.28
post_analysis_sparc0M
## [[1]]
## response0M_sparc ~ 0 + condition + (1 | condition:participant) + (1 | condition:trial)
##
## [[2]]
##
## [[3]]
## R2 SD CI CI_low CI_high CI_method Component Effectsize
## 1 0.82 0.02 0.95 0.78 0.85 HDI conditional Bayesian R-squared
## 2 0.10 0.11 0.95 0.00 0.34 HDI marginal Bayesian R-squared
##
## [[4]]
## Parameter Effects Component Mean CI CI_low
## 1 b_condition1 fixed conditional 10.90 0.95 10.67
## 2 b_condition2 fixed conditional 11.06 0.95 10.85
## 3 sd_condition:participant__Intercept random conditional 0.27 0.95 0.18
## 4 sd_condition:trial__Intercept random conditional 0.02 0.95 0.00
## 5 sigma fixed sigma 0.11 0.95 0.09
## CI_high pd log_BF Rhat ESS Group
## 1 11.06 1 77.94 1 615.23
## 2 11.23 1 58.30 1 432.37
## 3 0.45 1 16.69 1 404.60 condition:participant
## 4 0.07 1 -4.53 1 1290.88 condition:trial
## 5 0.14 1 33.15 1 1673.47
load(file = "Results/post_analysis_procustus1M.RData")
load(file = "Results/post_analysis_sparc1M.RData")
post_analysis_procustus1M
## [[1]]
## response1M_procustus ~ 0 + condition + (1 | condition:participant) + (1 | condition:trial)
##
## [[2]]
##
## [[3]]
## R2 SD CI CI_low CI_high CI_method Component Effectsize
## 1 0.51 0.09 0.95 0.32 0.64 HDI conditional Bayesian R-squared
## 2 0.06 0.07 0.95 0.00 0.24 HDI marginal Bayesian R-squared
##
## [[4]]
## Parameter Effects Component Mean CI CI_low
## 1 b_condition1 fixed conditional 0.17 0.95 0.04
## 2 b_condition2 fixed conditional 0.04 0.95 -0.10
## 3 sd_condition:participant__Intercept random conditional 0.16 0.95 0.09
## 4 sd_condition:trial__Intercept random conditional 0.08 0.95 0.02
## 5 sigma fixed sigma 0.15 0.95 0.13
## CI_high pd log_BF Rhat ESS Group
## 1 0.32 0.99 -0.17 1.00 851.47
## 2 0.18 0.71 -3.16 1.00 1022.25
## 3 0.24 1.00 3.55 1.01 609.55 condition:participant
## 4 0.19 1.00 -0.72 1.00 883.66 condition:trial
## 5 0.19 1.00 25.96 1.00 1336.37
post_analysis_sparc1M
## [[1]]
## response1M_sparc ~ 0 + condition + (1 | condition:participant) + (1 | condition:trial)
##
## [[2]]
##
## [[3]]
## R2 SD CI CI_low CI_high CI_method Component Effectsize
## 1 0.62 0.05 0.95 0.5 0.70 HDI conditional Bayesian R-squared
## 2 0.21 0.11 0.95 0.0 0.38 HDI marginal Bayesian R-squared
##
## [[4]]
## Parameter Effects Component Mean CI CI_low
## 1 b_condition1 fixed conditional -0.09 0.95 -0.17
## 2 b_condition2 fixed conditional 0.07 0.95 -0.01
## 3 sd_condition:participant__Intercept random conditional 0.13 0.95 0.08
## 4 sd_condition:trial__Intercept random conditional 0.02 0.95 0.00
## 5 sigma fixed sigma 0.11 0.95 0.09
## CI_high pd log_BF Rhat ESS Group
## 1 0.01 0.96 -1.93 1 735.74
## 2 0.16 0.96 -2.17 1 819.52
## 3 0.19 1.00 8.20 1 632.93 condition:participant
## 4 0.07 1.00 -3.57 1 1204.77 condition:trial
## 5 0.13 1.00 32.02 1 1934.52
We show the contrasts both for the procustus model_1 and the sparc
model_1. See the paper for a discussion about the contrast results. The
Label is coded as follows: c stands for condition, c12 for
a contrast of condition 1 and condition 2. t stands for trial, t12
stands for a contrast of trial 1 and trial 2. Accordingly, c1t12 stands
for a contrast of trial 1 an trial 2 in condition 1. c12t1 stands for a
contrast of trial 1 in condiition 1 versus condition 2.
load(file = "Results/hypothesis_test_procustus0M.RData")
load(file = "Results/hypothesis_test_sparc0M.RData")
hypothesis_test_procustus0M[[1]][,1]
## [1] "(b_condition1)-(b_condition2) < 0"
## [2] "(r_condition:trial[1_1,Intercept])-(r_condition:trial[1_2,Intercept]) > 0"
## [3] "(r_condition:trial[1_1,Intercept])-(r_condition:trial[1_3,Intercept]) > 0"
## [4] "(r_condition:trial[1_1,Intercept])-(r_condition:trial[1_4,Intercept]) > 0"
## [5] "(r_condition:trial[1_2,Intercept])-(r_condition:trial[1_3,Intercept]) > 0"
## [6] "(r_condition:trial[1_2,Intercept])-(r_condition:trial[1_4,Intercept]) > 0"
## [7] "(r_condition:trial[1_3,Intercept])-(r_condition:trial[1_4,Intercept]) > 0"
## [8] "(r_condition:trial[2_1,Intercept])-(r_condition:trial[2_2,Intercept]) > 0"
## [9] "(r_condition:trial[2_1,Intercept])-(r_condition:trial[2_3,Intercept]) > 0"
## [10] "(r_condition:trial[2_1,Intercept])-(r_condition:trial[2_4,Intercept]) > 0"
## [11] "(r_condition:trial[2_2,Intercept])-(r_condition:trial[2_3,Intercept]) > 0"
## [12] "(r_condition:trial[2_2,Intercept])-(r_condition:trial[2_4,Intercept]) > 0"
## [13] "(r_condition:trial[2_3,Intercept])-(r_condition:trial[2_4,Intercept]) > 0"
## [14] "(b_condition1+r_condition:trial[1_1,Intercept])-(b_condition2+r_condition:trial[2_1,Intercept]) > 0"
## [15] "(b_condition1+r_condition:trial[1_2,Intercept])-(b_condition2+r_condition:trial[2_2,Intercept]) > 0"
## [16] "(b_condition1+r_condition:trial[1_3,Intercept])-(b_condition2+r_condition:trial[2_3,Intercept]) > 0"
## [17] "(b_condition1+r_condition:trial[1_4,Intercept])-(b_condition2+r_condition:trial[2_4,Intercept]) > 0"
hypothesis_test_procustus0M[[1]][,-1]
## Label Estimate CI.Lower CI.Upper Post.Prob Star
## t1 c12 0.13 -0.02 0.28 0.07
## t2 c1t12 0.05 -0.03 0.14 0.83
## t3 c1t13 0.11 0.01 0.21 0.97 *
## t4 c1t14 0.14 0.04 0.24 0.99 *
## t5 c1t23 0.06 -0.02 0.14 0.90
## t6 c1t24 0.10 0.01 0.18 0.98 *
## t7 c1t34 0.04 -0.04 0.11 0.80
## t8 c2t12 0.06 -0.02 0.15 0.88
## t9 c2t13 0.06 -0.01 0.14 0.92
## t10 c2t14 0.07 0.00 0.14 0.93
## t11 c2t23 0.00 -0.08 0.08 0.54
## t12 c2t24 0.01 -0.07 0.09 0.59
## t13 c2t34 0.01 -0.07 0.08 0.55
## t14 c12t1 0.16 0.02 0.30 0.97 *
## t15 c12t2 0.17 0.03 0.32 0.97 *
## t16 c12t3 0.12 -0.02 0.26 0.92
## t17 c12t4 0.09 -0.06 0.22 0.84
hypothesis_test_procustus0M[[2]]
hypothesis_test_sparc0M[[1]][,1]
## [1] "(b_condition1)-(b_condition2) < 0"
## [2] "(r_condition:trial[1_1,Intercept])-(r_condition:trial[1_2,Intercept]) > 0"
## [3] "(r_condition:trial[1_1,Intercept])-(r_condition:trial[1_3,Intercept]) > 0"
## [4] "(r_condition:trial[1_1,Intercept])-(r_condition:trial[1_4,Intercept]) > 0"
## [5] "(r_condition:trial[1_2,Intercept])-(r_condition:trial[1_3,Intercept]) > 0"
## [6] "(r_condition:trial[1_2,Intercept])-(r_condition:trial[1_4,Intercept]) > 0"
## [7] "(r_condition:trial[1_3,Intercept])-(r_condition:trial[1_4,Intercept]) > 0"
## [8] "(r_condition:trial[2_1,Intercept])-(r_condition:trial[2_2,Intercept]) > 0"
## [9] "(r_condition:trial[2_1,Intercept])-(r_condition:trial[2_3,Intercept]) > 0"
## [10] "(r_condition:trial[2_1,Intercept])-(r_condition:trial[2_4,Intercept]) > 0"
## [11] "(r_condition:trial[2_2,Intercept])-(r_condition:trial[2_3,Intercept]) > 0"
## [12] "(r_condition:trial[2_2,Intercept])-(r_condition:trial[2_4,Intercept]) > 0"
## [13] "(r_condition:trial[2_3,Intercept])-(r_condition:trial[2_4,Intercept]) > 0"
## [14] "(b_condition1+r_condition:trial[1_1,Intercept])-(b_condition2+r_condition:trial[2_1,Intercept]) > 0"
## [15] "(b_condition1+r_condition:trial[1_2,Intercept])-(b_condition2+r_condition:trial[2_2,Intercept]) > 0"
## [16] "(b_condition1+r_condition:trial[1_3,Intercept])-(b_condition2+r_condition:trial[2_3,Intercept]) > 0"
## [17] "(b_condition1+r_condition:trial[1_4,Intercept])-(b_condition2+r_condition:trial[2_4,Intercept]) > 0"
hypothesis_test_sparc0M[[1]][,-1]
## Label Estimate CI.Lower CI.Upper Post.Prob Star
## t1 c12 -0.16 -0.36 0.03 0.92
## t2 c1t12 0.01 -0.02 0.06 0.65
## t3 c1t13 0.01 -0.03 0.05 0.61
## t4 c1t14 0.00 -0.03 0.05 0.55
## t5 c1t23 0.00 -0.05 0.03 0.45
## t6 c1t24 -0.01 -0.05 0.03 0.39
## t7 c1t34 0.00 -0.04 0.03 0.45
## t8 c2t12 -0.01 -0.06 0.02 0.32
## t9 c2t13 -0.01 -0.05 0.03 0.40
## t10 c2t14 0.01 -0.02 0.06 0.65
## t11 c2t23 0.01 -0.03 0.05 0.59
## t12 c2t24 0.02 -0.01 0.09 0.79
## t13 c2t34 0.02 -0.02 0.08 0.73
## t14 c12t1 -0.16 -0.35 0.04 0.09
## t15 c12t2 -0.18 -0.38 0.01 0.06
## t16 c12t3 -0.17 -0.37 0.02 0.07
## t17 c12t4 -0.15 -0.35 0.04 0.10
hypothesis_test_sparc0M[[2]]
load(file = "Results/hypothesis_test_procustus1M.RData")
load(file = "Results/hypothesis_test_sparc1M.RData")
hypothesis_test_procustus1M[[1]][,1]
## [1] "(b_condition1)-(b_condition2) < 0"
## [2] "(r_condition:trial[1_1,Intercept])-(r_condition:trial[1_2,Intercept]) > 0"
## [3] "(r_condition:trial[1_1,Intercept])-(r_condition:trial[1_3,Intercept]) > 0"
## [4] "(r_condition:trial[1_1,Intercept])-(r_condition:trial[1_4,Intercept]) > 0"
## [5] "(r_condition:trial[1_2,Intercept])-(r_condition:trial[1_3,Intercept]) > 0"
## [6] "(r_condition:trial[1_2,Intercept])-(r_condition:trial[1_4,Intercept]) > 0"
## [7] "(r_condition:trial[1_3,Intercept])-(r_condition:trial[1_4,Intercept]) > 0"
## [8] "(r_condition:trial[2_1,Intercept])-(r_condition:trial[2_2,Intercept]) > 0"
## [9] "(r_condition:trial[2_1,Intercept])-(r_condition:trial[2_3,Intercept]) > 0"
## [10] "(r_condition:trial[2_1,Intercept])-(r_condition:trial[2_4,Intercept]) > 0"
## [11] "(r_condition:trial[2_2,Intercept])-(r_condition:trial[2_3,Intercept]) > 0"
## [12] "(r_condition:trial[2_2,Intercept])-(r_condition:trial[2_4,Intercept]) > 0"
## [13] "(r_condition:trial[2_3,Intercept])-(r_condition:trial[2_4,Intercept]) > 0"
## [14] "(b_condition1+r_condition:trial[1_1,Intercept])-(b_condition2+r_condition:trial[2_1,Intercept]) > 0"
## [15] "(b_condition1+r_condition:trial[1_2,Intercept])-(b_condition2+r_condition:trial[2_2,Intercept]) > 0"
## [16] "(b_condition1+r_condition:trial[1_3,Intercept])-(b_condition2+r_condition:trial[2_3,Intercept]) > 0"
## [17] "(b_condition1+r_condition:trial[1_4,Intercept])-(b_condition2+r_condition:trial[2_4,Intercept]) > 0"
hypothesis_test_procustus1M[[1]][,-1]
## Label Estimate CI.Lower CI.Upper Post.Prob Star
## t1 c12 0.13 -0.02 0.29 0.08
## t2 c1t12 0.06 -0.02 0.15 0.86
## t3 c1t13 0.12 0.02 0.21 0.98 *
## t4 c1t14 0.16 0.05 0.26 0.99 *
## t5 c1t23 0.06 -0.02 0.14 0.89
## t6 c1t24 0.10 0.01 0.18 0.98 *
## t7 c1t34 0.04 -0.03 0.11 0.80
## t8 c2t12 0.06 -0.01 0.15 0.90
## t9 c2t13 0.07 -0.01 0.15 0.92
## t10 c2t14 0.07 0.00 0.15 0.94
## t11 c2t23 0.00 -0.08 0.08 0.55
## t12 c2t24 0.00 -0.08 0.08 0.55
## t13 c2t34 0.00 -0.08 0.08 0.52
## t14 c12t1 0.16 0.02 0.31 0.97 *
## t15 c12t2 0.17 0.03 0.31 0.98 *
## t16 c12t3 0.11 -0.02 0.25 0.91
## t17 c12t4 0.08 -0.06 0.21 0.83
hypothesis_test_procustus1M[[2]]
hypothesis_test_sparc1M[[1]][,1]
## [1] "(b_condition1)-(b_condition2) < 0"
## [2] "(r_condition:trial[1_1,Intercept])-(r_condition:trial[1_2,Intercept]) > 0"
## [3] "(r_condition:trial[1_1,Intercept])-(r_condition:trial[1_3,Intercept]) > 0"
## [4] "(r_condition:trial[1_1,Intercept])-(r_condition:trial[1_4,Intercept]) > 0"
## [5] "(r_condition:trial[1_2,Intercept])-(r_condition:trial[1_3,Intercept]) > 0"
## [6] "(r_condition:trial[1_2,Intercept])-(r_condition:trial[1_4,Intercept]) > 0"
## [7] "(r_condition:trial[1_3,Intercept])-(r_condition:trial[1_4,Intercept]) > 0"
## [8] "(r_condition:trial[2_1,Intercept])-(r_condition:trial[2_2,Intercept]) > 0"
## [9] "(r_condition:trial[2_1,Intercept])-(r_condition:trial[2_3,Intercept]) > 0"
## [10] "(r_condition:trial[2_1,Intercept])-(r_condition:trial[2_4,Intercept]) > 0"
## [11] "(r_condition:trial[2_2,Intercept])-(r_condition:trial[2_3,Intercept]) > 0"
## [12] "(r_condition:trial[2_2,Intercept])-(r_condition:trial[2_4,Intercept]) > 0"
## [13] "(r_condition:trial[2_3,Intercept])-(r_condition:trial[2_4,Intercept]) > 0"
## [14] "(b_condition1+r_condition:trial[1_1,Intercept])-(b_condition2+r_condition:trial[2_1,Intercept]) > 0"
## [15] "(b_condition1+r_condition:trial[1_2,Intercept])-(b_condition2+r_condition:trial[2_2,Intercept]) > 0"
## [16] "(b_condition1+r_condition:trial[1_3,Intercept])-(b_condition2+r_condition:trial[2_3,Intercept]) > 0"
## [17] "(b_condition1+r_condition:trial[1_4,Intercept])-(b_condition2+r_condition:trial[2_4,Intercept]) > 0"
hypothesis_test_sparc1M[[1]][,-1]
## Label Estimate CI.Lower CI.Upper Post.Prob Star
## t1 c12 -0.16 -0.25 -0.06 0.99 *
## t2 c1t12 0.01 -0.02 0.07 0.69
## t3 c1t13 0.01 -0.03 0.05 0.60
## t4 c1t14 0.00 -0.03 0.05 0.57
## t5 c1t23 -0.01 -0.05 0.03 0.41
## t6 c1t24 -0.01 -0.06 0.02 0.37
## t7 c1t34 0.00 -0.04 0.03 0.45
## t8 c2t12 -0.01 -0.05 0.02 0.34
## t9 c2t13 -0.01 -0.05 0.03 0.42
## t10 c2t14 0.02 -0.02 0.07 0.70
## t11 c2t23 0.00 -0.03 0.05 0.57
## t12 c2t24 0.03 -0.01 0.09 0.79
## t13 c2t34 0.02 -0.02 0.08 0.75
## t14 c12t1 -0.15 -0.25 -0.05 0.01
## t15 c12t2 -0.18 -0.28 -0.07 0.00
## t16 c12t3 -0.17 -0.27 -0.06 0.01
## t17 c12t4 -0.14 -0.24 -0.04 0.02
hypothesis_test_sparc1M[[2]]
In Figure 5 and 6 we show the trial distributions drawn from the
calibrated model for the procustus and sparc metrics. The graphs below
are based on the calibrated models from which we took the group-level
effects of trial and build their posterior distributions.
In this approach, as known, estimates are typically shrinked towards the
mean, assuring robust modelling. Alternatively, we could have modelled a
longitudinal model using time (weeks or days) as temporal variable.
However, in this context, we believe that order was more relevant than
time and therefore we coded trial as a factor rather than a numeric
variable.
Part 2 of this analysis is related to questionnaire models and hypothesis 3 and 4, following the workflow of Figure 1 (bottom part).
We tested the models for the 4 questions and report here the log of the Bayes factor.
load(file = "Results/BF_model_comparison_WPQ_0M.RData")
load(file = "Results/BF_model_comparison_MPQS_0M.RData")
load(file = "Results/BF_model_comparison_MPQP_0M.RData")
load(file = "Results/BF_model_comparison_Difficulty_0M.RData")
print(paste("Bayes factor in favor of model_3 over model_4 (WPQ_0M): ", BF_model_comparison_WPQ_0M ))
## [1] "Bayes factor in favor of model_3 over model_4 (WPQ_0M): 3.07659020662719"
## [2] "Bayes factor in favor of model_3 over model_4 (WPQ_0M): TRUE"
print(paste("Bayes factor in favor of model_3 over model_4 (MPQS_0M): ", BF_model_comparison_MPQS_0M ))
## [1] "Bayes factor in favor of model_3 over model_4 (MPQS_0M): 3.42336468424023"
## [2] "Bayes factor in favor of model_3 over model_4 (MPQS_0M): TRUE"
print(paste("Bayes factor in favor of model_3 over model_4 (MPQP_0M): ", BF_model_comparison_MPQP_0M ))
## [1] "Bayes factor in favor of model_3 over model_4 (MPQP_0M): 3.27297564793447"
## [2] "Bayes factor in favor of model_3 over model_4 (MPQP_0M): TRUE"
print(paste("Bayes factor in favor of model_3 over model_4 (Difficulty_OM): ", BF_model_comparison_Difficulty_0M ))
## [1] "Bayes factor in favor of model_3 over model_4 (Difficulty_OM): 2.48310653516731"
## [2] "Bayes factor in favor of model_3 over model_4 (Difficulty_OM): TRUE"
load(file = "Results/BF_model_comparison_WPQ_1M.RData")
load(file = "Results/BF_model_comparison_MPQS_1M.RData")
load(file = "Results/BF_model_comparison_MPQP_1M.RData")
load(file = "Results/BF_model_comparison_Difficulty_1M.RData")
print(paste("Bayes factor in favor of model_3 over model_4 (WPQ_1M): ", BF_model_comparison_WPQ_1M ))
## [1] "Bayes factor in favor of model_3 over model_4 (WPQ_1M): 4.11902706997761"
## [2] "Bayes factor in favor of model_3 over model_4 (WPQ_1M): TRUE"
print(paste("Bayes factor in favor of model_3 over model_4 (MPQS_1M): ", BF_model_comparison_MPQS_1M ))
## [1] "Bayes factor in favor of model_3 over model_4 (MPQS_1M): 2.81496110400582"
## [2] "Bayes factor in favor of model_3 over model_4 (MPQS_1M): TRUE"
print(paste("Bayes factor in favor of model_3 over model_4 (MPQP_1M): ", BF_model_comparison_MPQP_1M ))
## [1] "Bayes factor in favor of model_3 over model_4 (MPQP_1M): 3.45789108826106"
## [2] "Bayes factor in favor of model_3 over model_4 (MPQP_1M): TRUE"
print(paste("Bayes factor in favor of model_3 over model_4 (Difficulty_1M): ", BF_model_comparison_Difficulty_1M ))
## [1] "Bayes factor in favor of model_3 over model_4 (Difficulty_1M): -0.241128804241981"
## [2] "Bayes factor in favor of model_3 over model_4 (Difficulty_1M): TRUE"
From this analysis it can be concluded that there is moderate
evidence for model_3 for presence, and anecdotical evidence for model_3
for difficulty. Basically, it means that procustus does not
contribute to an explanation of those responses.
We show the diagnostics both for the 4 question models:
the model_1 formula,
the Bayes_R2 analysis,
the model parameters,
the posterior prediction check (pp_check) next to the plot of
fixed parameters (i.e. condition)
Given the fact that procustus has no substantial
influence we use model_1 (without procustus)
load(file = "Results/post_analysis_procustus_basicWPQ.RData")
load(file = "Results/post_analysis_procustus_basicMPQS.RData")
load(file = "Results/post_analysis_procustus_basicMPQP.RData")
load(file = "Results/post_analysis_procustus_basicDifficulty.RData")
post_analysis_procustus_basicWPQ
## [[1]]
## WPQ ~ 0 + condition + (1 | condition:participant + condition:trial)
##
## [[2]]
##
## [[3]]
## R2 SD CI CI_low CI_high CI_method Component Effectsize
## 1 0.73 0.03 0.95 0.66 0.80 HDI conditional Bayesian R-squared
## 2 0.12 0.12 0.95 0.00 0.33 HDI marginal Bayesian R-squared
##
## [[4]]
## Parameter Effects Component Mean CI CI_low
## 1 b_condition1 fixed conditional 4.69 0.95 4.35
## 2 b_condition2 fixed conditional 5.06 0.95 4.74
## 3 sd_condition:participant__Intercept random conditional 0.49 0.95 0.34
## 4 sd_condition:trial__Intercept random conditional 0.08 0.95 0.00
## 5 sigma fixed sigma 0.30 0.95 0.25
## CI_high pd log_BF Rhat ESS Group
## 1 5.03 1 49.68 1 582.83
## 2 5.38 1 63.15 1 612.48
## 3 0.69 1 13.99 1 609.76 condition:participant
## 4 0.23 1 -3.02 1 1214.38 condition:trial
## 5 0.36 1 20.72 1 1823.86
post_analysis_procustus_basicMPQS
## [[1]]
## MPQS ~ 0 + condition + (1 | condition:participant + condition:trial)
##
## [[2]]
##
## [[3]]
## R2 SD CI CI_low CI_high CI_method Component Effectsize
## 1 0.66 0.05 0.95 0.56 0.74 HDI conditional Bayesian R-squared
## 2 0.17 0.11 0.95 0.00 0.36 HDI marginal Bayesian R-squared
##
## [[4]]
## Parameter Effects Component Mean CI CI_low
## 1 b_condition1 fixed conditional 2.74 0.95 2.25
## 2 b_condition2 fixed conditional 3.55 0.95 3.08
## 3 sd_condition:participant__Intercept random conditional 0.72 0.95 0.49
## 4 sd_condition:trial__Intercept random conditional 0.16 0.95 0.01
## 5 sigma fixed sigma 0.58 0.95 0.48
## CI_high pd log_BF Rhat ESS Group
## 1 3.25 1 17.97 1.01 821.75
## 2 4.04 1 24.57 1.00 909.62
## 3 1.05 1 10.69 1.00 758.32 condition:participant
## 4 0.45 1 -1.34 1.00 908.25 condition:trial
## 5 0.69 1 32.41 1.00 2631.10
post_analysis_procustus_basicMPQP
## [[1]]
## MPQP ~ 0 + condition + (1 | condition:participant + condition:trial)
##
## [[2]]
##
## [[3]]
## R2 SD CI CI_low CI_high CI_method Component Effectsize
## 1 0.65 0.05 0.95 0.54 0.74 HDI conditional Bayesian R-squared
## 2 0.23 0.13 0.95 0.00 0.41 HDI marginal Bayesian R-squared
##
## [[4]]
## Parameter Effects Component Mean CI CI_low
## 1 b_condition1 fixed conditional 3.49 0.95 2.99
## 2 b_condition2 fixed conditional 4.42 0.95 3.88
## 3 sd_condition:participant__Intercept random conditional 0.68 0.95 0.45
## 4 sd_condition:trial__Intercept random conditional 0.20 0.95 0.01
## 5 sigma fixed sigma 0.59 0.95 0.49
## CI_high pd log_BF Rhat ESS Group
## 1 3.98 1 29.71 1 934.40
## 2 4.92 1 34.57 1 708.05
## 3 1.00 1 9.99 1 690.54 condition:participant
## 4 0.52 1 -0.94 1 1120.26 condition:trial
## 5 0.71 1 35.85 1 1605.05
post_analysis_procustus_basicDifficulty
## [[1]]
## Difficulty ~ 0 + condition + (1 | condition:participant + condition:trial)
##
## [[2]]
##
## [[3]]
## R2 SD CI CI_low CI_high CI_method Component Effectsize
## 1 0.71 0.04 0.95 0.62 0.77 HDI conditional Bayesian R-squared
## 2 0.03 0.03 0.95 0.00 0.15 HDI marginal Bayesian R-squared
##
## [[4]]
## Parameter Effects Component Mean CI CI_low
## 1 b_condition1 fixed conditional 0.10 0.95 -0.49
## 2 b_condition2 fixed conditional -0.11 0.95 -0.66
## 3 sd_condition:participant__Intercept random conditional 0.89 0.95 0.62
## 4 sd_condition:trial__Intercept random conditional 0.11 0.95 0.00
## 5 sigma fixed sigma 0.54 0.95 0.46
## CI_high pd log_BF Rhat ESS Group
## 1 0.64 0.66 -2.33 1 421.17
## 2 0.44 0.66 -2.27 1 458.74
## 3 1.29 1.00 13.90 1 660.84 condition:participant
## 4 0.34 1.00 -2.69 1 1422.96 condition:trial
## 5 0.66 1.00 25.24 1 1824.16
load(file = "Results/Contrast_WPQ.RData")
load(file = "Results/Contrast_MPQS.RData")
load(file = "Results/Contrast_MPQP.RData")
load(file = "Results/Contrast_Difficulty.RData")
Contrasts_Questionnaires <- rbind(
Contrast_WPQ[[1]][-1],
Contrast_MPQS[[1]][-1],
Contrast_WPQP[[1]][-1],
Contrast_Difficulty[[1]][-1]
)
rownames(Contrasts_Questionnaires) <- NULL
Contrasts_Questionnaires
## Label Estimate CI.Lower CI.Upper Post.Prob Star
## 1 c12 -0.37 -0.77 0.02 0.94
## 2 c1t12 0.05 -0.08 0.21 0.70
## 3 c1t13 0.06 -0.06 0.23 0.75
## 4 c1t14 0.02 -0.10 0.17 0.61
## 5 c1t23 0.02 -0.11 0.16 0.57
## 6 c1t24 -0.02 -0.18 0.11 0.39
## 7 c1t34 -0.04 -0.19 0.09 0.34
## 8 c2t12 -0.10 -0.28 0.03 0.17
## 9 c2t13 -0.04 -0.20 0.08 0.32
## 10 c2t14 -0.12 -0.33 0.02 0.13
## 11 c2t23 0.05 -0.07 0.22 0.71
## 12 c2t24 -0.02 -0.16 0.11 0.41
## 13 c2t34 -0.07 -0.25 0.05 0.22
## 14 c12t1 -0.28 -0.67 0.11 0.12
## 15 c12t2 -0.42 -0.82 -0.03 0.04
## 16 c12t3 -0.38 -0.77 0.01 0.05
## 17 c12t4 -0.42 -0.82 -0.03 0.04
## 18 c12 -0.80 -1.38 -0.24 0.98 *
## 19 c1t12 0.17 -0.08 0.56 0.80
## 20 c1t13 0.23 -0.03 0.65 0.87
## 21 c1t14 0.17 -0.07 0.53 0.81
## 22 c1t23 0.06 -0.18 0.36 0.64
## 23 c1t24 0.00 -0.27 0.26 0.49
## 24 c1t34 -0.07 -0.36 0.17 0.35
## 25 c2t12 -0.07 -0.38 0.17 0.37
## 26 c2t13 0.03 -0.23 0.30 0.58
## 27 c2t14 0.08 -0.16 0.39 0.68
## 28 c2t23 0.10 -0.13 0.41 0.71
## 29 c2t24 0.15 -0.08 0.50 0.78
## 30 c2t34 0.05 -0.18 0.32 0.62
## 31 c12t1 -0.67 -1.25 -0.05 0.04
## 32 c12t2 -0.91 -1.50 -0.33 0.01
## 33 c12t3 -0.87 -1.46 -0.29 0.01
## 34 c12t4 -0.75 -1.31 -0.17 0.02
## 35 c12 -0.93 -1.51 -0.36 0.99 *
## 36 c1t12 -0.03 -0.35 0.26 0.45
## 37 c1t13 0.23 -0.06 0.62 0.86
## 38 c1t14 0.04 -0.25 0.34 0.59
## 39 c1t23 0.26 -0.05 0.68 0.88
## 40 c1t24 0.07 -0.22 0.39 0.63
## 41 c1t34 -0.19 -0.56 0.09 0.17
## 42 c2t12 -0.21 -0.60 0.08 0.16
## 43 c2t13 0.04 -0.24 0.33 0.58
## 44 c2t14 -0.04 -0.36 0.25 0.40
## 45 c2t23 0.25 -0.06 0.66 0.87
## 46 c2t24 0.17 -0.10 0.53 0.80
## 47 c2t34 -0.08 -0.41 0.19 0.33
## 48 c12t1 -0.82 -1.39 -0.24 0.01
## 49 c12t2 -1.00 -1.56 -0.42 0.00
## 50 c12t3 -1.01 -1.60 -0.43 0.00
## 51 c12t4 -0.90 -1.45 -0.33 0.01
## 52 c12 0.22 -0.46 0.84 0.28
## 53 c1t12 0.03 -0.15 0.26 0.59
## 54 c1t13 0.09 -0.08 0.38 0.72
## 55 c1t14 0.00 -0.21 0.21 0.49
## 56 c1t23 0.06 -0.11 0.32 0.66
## 57 c1t24 -0.03 -0.27 0.16 0.41
## 58 c1t34 -0.09 -0.40 0.09 0.27
## 59 c2t12 0.05 -0.12 0.30 0.66
## 60 c2t13 0.02 -0.18 0.24 0.56
## 61 c2t14 -0.02 -0.24 0.17 0.44
## 62 c2t23 -0.03 -0.26 0.15 0.40
## 63 c2t24 -0.07 -0.35 0.10 0.31
## 64 c2t34 -0.04 -0.28 0.15 0.39
## 65 c12t1 0.24 -0.42 0.88 0.72
## 66 c12t2 0.26 -0.39 0.90 0.75
## 67 c12t3 0.17 -0.50 0.80 0.67
## 68 c12t4 0.22 -0.45 0.87 0.71